# 导入需要的模块
import pandas as pd
import matplotlib.pyplot as plt
import numpy as np
# 用来绘图的,封装了matplot
# 要注意的是一旦导入了seaborn,
# matplotlib的默认作图风格就会被覆盖成seaborn的格式
import seaborn as sns
from scipy import stats
from scipy.stats import norm
from sklearn.preprocessing import StandardScaler
import warnings
warnings.filterwarnings('ignore')
%matplotlib inline
# 为了在jupyter notebook里作图,需要用到这个命令
from mxnet import autograd, gluon, init, nd
from mxnet.gluon import data as gdata, loss as gloss, nn
from IPython import display
data_train = pd.read_csv("./data/train.csv")
# data_train
因为我们最后关心的房价,首先直接来看一下房价分布情况:
绘制出关于房价的直方图和密度曲线,可见房价的密度分布接近正态分布,进一步可以使用峰度(Kurtosis)和 偏度(Skewness)来描述。
峰度(Kurtosis):描述某变量所有取值分布形态陡缓程度的统计量。 它是和正态分布相比较的。
Kurtosis=0 与正态分布的陡缓程度相同。
Kurtosis>0 比正态分布的高峰更加陡峭——尖顶峰
Kurtosis<0 比正态分布的高峰来得平坦——平顶峰计算公式:β = M_4 /σ^4 偏度:
偏度(Skewness):描述某变量取值分布对称性的统计量。
Skewness=0 分布形态与正态分布偏度相同
Skewness>0 正偏差数值较大,为正偏或右偏。长尾巴拖在右边。
Skewness<0 负偏差数值较大,为负偏或左偏。长尾巴拖在左边。 计算公式: S= (X^ - M_0)/δ Skewness 越大,分布形态偏移程度越大。
可以发现:峰度和坡度均大于0,说明房价的概率密度曲线相比于正态分布高峰更加陡峭,向右边偏移;房价更多的集中在150000附近,并且150000的右侧概率相比于左侧概率更密集。
# data_train['SalePrice'].describe()
sns.distplot(data_train['SalePrice'])
# skewness and kurtosis
print(f"Skewness: {data_train['SalePrice'].skew()}")
print(f"Kurtosis: {data_train['SalePrice'].kurt()}")
然后我们来看看可能决定房价的其他属性,数据集中给出的有79个,将这79个信息进行分类、概括,以便进行特征的选取。
首先将属性类型分类为三种:“building--0”、“space--1”、“location--2”。具体解释如下:
然后每种属性的数据类型分类两类:“数值型--0”和“类别型--1”
分类后的结果如下:
| 变量名 | 中文解释 | 数据类型 | 分类 |
|---|---|---|---|
| MSSubClass | 住宅类型 | 1 | building |
| MSZoning | 住宅所在分区类型 | 1 | location |
| LotFrontage | 房子到街区的距离 | 0 | location |
| LotArea | 面积 | 0 | space |
| Street | 所在街道 | 1 | location |
| Alley | 所在小路 | 1 | location |
| LotShape | 户型不规则程度 | 1 | building |
| LandContour | 平坦程度 | 1 | building |
| Utilities | 水电气供应 | 1 | building |
| LotConfig | 位置(临街状态) | 1 | location |
| LandSlope | 倾斜程度 | 1 | building |
| Neighborhood | 邻接情况 | 1 | location |
| Condition1 | 位置特征1(道路、铁路) | 1 | location |
| Condition2 | 位置特征2(道路、铁路) | 1 | location |
| BldgType | 住宅类型(楼房类型) | 1 | building |
| HouseStyle | 风格 | 1 | building |
| OverallQual | 房屋的整体情况评级 | 1 | building |
| OverallCond | 房屋的整体状况评级 | 1 | building |
| YearBuilt | 建造年份 | 1 | building |
| YearRemodAdd | 重修年份 | 1 | building |
| RoofStyle | 屋顶类型 | 1 | building |
| RoofMatl | 屋顶材料 | 1 | building |
| Exterior1st | 外墙材料1 | 1 | building |
| Exterior2nd | 外墙材料2 | 1 | building |
| MasVnrType | 石材类型 | 1 | building |
| MasVnrArea | 石材面积 | 0 | space |
| ExterQual | 外部材料质量评级 | 1 | building |
| ExterCond | 外部材料现状评级 | 1 | building |
| Foundation | 地基类型 | 1 | building |
| BsmtQual | 地基深度 | 0 | space |
| BsmtCond | 地基现状评级 | 1 | building |
| BsmtExposure | 地下室出口暴露情况 | 1 | building |
| BsmtFinType1 | 地下室第一次状况评级 | 1 | building |
| BsmtFinSF1 | 地下室第一次完工面积 | 0 | space |
| BsmtFinType2 | 地下室第二次状况评级 | 1 | building |
| BsmtFinSF2 | 地下室第二次完工面积 | 0 | space |
| BsmtUnfSF | 地下室未完成面积 | 0 | space |
| TotalBsmtSF | 地下室总面积 | 0 | space |
| Heating | 供暖类型 | 1 | building |
| HeatingQC | 供热质量状况评级 | 1 | building |
| CentralAir | 是否有中央空调 | 1 | building |
| Electrical | 电力系统类型(熔断机制) | 1 | building |
| 1stFlrSF | 一楼面积 | 0 | space |
| 2ndFlrSF | 二楼面积 | 0 | space |
| LowQualFinSF | 低质量完工面积 | 0 | space |
| GrLivArea | 地面以上面积 | 0 | space |
| BsmtFullBath | 地下室全浴室 | 0 | space |
| BsmtHalfBath | 地下室半浴室 | 0 | space |
| FullBath | 地面以上全浴室 | 0 | space |
| HalfBath | 地面以上半浴室 | 0 | space |
| BedroomAbvGr | 地面以上卧室 | 0 | space |
| KitchenAbvGr | 地面以上厨房 | 0 | space |
| KitchenQual | 厨房质量评级 | 1 | building |
| TotRmsAbvGrd | 地面以上房间数 | 0 | space |
| Functional | 房屋功能完整程度 | 1 | building |
| Fireplaces | 壁炉数量 | 1 | building |
| FireplaceQu | 壁炉质量评级 | 1 | building |
| GarageType | 车库类型 | 1 | building |
| GarageYrBlt | 车库建造年份 | 1 | building |
| GarageFinish | 车库装修是否完成 | 1 | building |
| GarageCars | 车库容量 | 0 | space |
| GarageArea | 车库面积 | 0 | space |
| GarageQual | 车库质量评级 | 1 | building |
| GarageCond | 车库现状评级 | 1 | building |
| PavedDrive | 车道铺装类型 | 1 | building |
| WoodDeckSF | 木质地板面积 | 0 | space |
| OpenPorchSF | 开放式门廊面积 | 0 | space |
| EnclosedPorch | 封闭式门廊面积 | 0 | space |
| 3SsnPorch | 时令门廊面积 | 0 | space |
| ScreenPorch | 屏幕门廊面积 | 0 | space |
| PoolArea | 游泳池面积 | 0 | space |
| PoolQC | 游泳池质量评级 | 1 | building |
| Fence | 栅栏质量评级 | 1 | building |
| MiscFeature | 其他功能 | 1 | building |
| MiscVal | 该其他功能价值 | 1 | building |
| MoSold | 卖出月份 | 1 | building |
| YrSold | 卖出年份 | 1 | building |
| SaleType | 交易类型 | 1 | building |
| SaleCondition | 交易情况 | 1 | building |
然后继续在这些特征中找一些认为比较重要的特征详细的看看。
有中央空调的房子售价普遍要高一些
# CentralAir
var = 'CentralAir'
data = pd.concat([data_train['SalePrice'], data_train[var]], axis=1) # 横向表的拼接
fig = sns.boxplot(x=var, y="SalePrice", data=data) # 绘制箱型图
fig.axis(ymin=0, ymax=800000);
评分越高售价越高
# OverallQual
var = 'OverallQual'
data = pd.concat([data_train['SalePrice'], data_train[var]], axis=1)
fig = sns.boxplot(x=var, y="SalePrice", data=data)
fig.axis(ymin=0, ymax=800000);
规律似乎不太明显,但是可以看出新建的房屋的售价会相对高一点,但是影响比较小。这一点散点图比箱型图表现的更加明显。
# YearBuilt boxplot
var = 'YearBuilt'
data = pd.concat([data_train['SalePrice'], data_train[var]], axis=1)
f, ax = plt.subplots(figsize=(26, 12)) # 指定图片大小
fig = sns.boxplot(x=var, y="SalePrice", data=data)
fig.axis(ymin=0, ymax=800000);
# YearBuilt scatter
var = 'YearBuilt'
data = pd.concat([data_train['SalePrice'], data_train[var]], axis=1)
data.plot.scatter(x=var, y="SalePrice", ylim=(0, 800000))
在一些地段(如NoRidge、NridgHt、StoneBr)房屋i的售价偏高,而另一些地段(如BrkSide、MeadowV)房屋的售价偏低;可见可对地段做一个聚类,存在高价格地段、低价格地段的分类。
# Neighborhood boxplot
var = 'Neighborhood'
data = pd.concat([data_train['SalePrice'], data_train[var]], axis=1)
f, ax = plt.subplots(figsize=(26, 12))
fig = sns.boxplot(x=var, y="SalePrice", data=data)
fig.axis(ymin=0, ymax=800000);
var ="LotArea"
data=pd.concat([data_train['SalePrice'],data_train[var]],axis=1)
data.plot.scatter(x=var,y='SalePrice',ylim=(0.800000))
地表面积的增加可能会导致房屋售价的增加,但是这种影响是不确定的,仍存在相当一部分数据,大地表面积的房屋和小地表面积的房屋售价相近。
但可以确认的是,地表面积影响了房屋售价的上限,房屋售价的上限大致与地表面积呈正比。
var ='GrLivArea'
data=pd.concat([data_train['SalePrice'],data_train[var]],axis=1)
data.plot.scatter(x=var,y='SalePrice',ylim=(0,800000))
地下室面积的增大会导致售价的升高,但这种影响不是一定的。
var ='TotalBsmtSF'
data=pd.concat([data_train['SalePrice'],data_train[var]],axis=1)
data.plot.scatter(x=var,y='SalePrice',ylim=(0,800000))
var='MiscVal'
data=pd.concat([data_train['SalePrice'],data_train[var]],axis=1)
data.plot.scatter(x=var,y='SalePrice',ylim=(0,800000))
车库面积的增大可能会导致售价增高
var='GarageArea'
data=pd.concat([data_train['SalePrice'],data_train[var]],axis=1)
data.plot.scatter(x=var,y='SalePrice',ylim=(0,800000))
车库容量的增大会导致房屋售价的上限增大。
var='GarageCars'
data=pd.concat([data_train['SalePrice'],data_train[var]],axis=1)
data.plot.scatter(x=var,y='SalePrice',ylim=(0,800000))
上面选择了一些自己关心的属性进行分析,但是这是不全面和相当主观的,下面会用关系矩阵和关系点图进行进一步分析。
在开始之前可以先用散点图的方式将数据的分布情况画出,让我们有一个基本的感受。
# 数值量特征
feats_numeric = data_train.dtypes[data_train.dtypes != "object"].index.values
#feats_numeric = [attr for attr in df_allX.columns if df_allX.dtypes[attr] != 'object']
# 字符量特征
feats_object = data_train.dtypes[data_train.dtypes == "object"].index.values
#feats_object = [attr for attr in df_allX.columns if df_allX.dtypes[attr] == 'object']
#feats_object = df_train.select_dtypes(include = ["object"]).columns
# 离散的数值量,需要人工甄别
feats_numeric_discrete = ['MSSubClass','OverallQual','OverallCond'] # 户型、整体质量打分、整体条件打分 —— 文档中明确定义的类型量
feats_numeric_discrete += ['TotRmsAbvGrd','KitchenAbvGr','BedroomAbvGr','GarageCars','Fireplaces'] # 房间数量
feats_numeric_discrete += ['FullBath','HalfBath','BsmtHalfBath','BsmtFullBath'] # 外国人这么爱洗澡?搞这么多浴室
feats_numeric_discrete += ['MoSold','YrSold'] # 年、月,这些不看成离散的应该也行
# 连续型特征
feats_continu = feats_numeric.copy()
# 离散型特征
feats_discrete = feats_object.copy()
for f in feats_numeric_discrete:
feats_continu = np.delete(feats_continu,np.where(feats_continu == f))
feats_discrete = np.append(feats_discrete,f)
print(feats_continu.shape,feats_discrete.shape)
71个特征中有26个来连续型特征,55个离散性特征。
def plotfeats(frame,feats,kind,cols=4):
"""批量绘图函数。
Parameters
----------
frame : pandas.DataFrame
待绘图的数据
feats : list 或 numpy.array
待绘图的列名称
kind : str
绘图格式:'hist'-直方图;'scatter'-散点图;'hs'-直方图和散点图隔行交替;'box'-箱线图,每个feat一幅图;'boxp'-Price做纵轴,feat做横轴的箱线图。
cols : int
每行绘制几幅图
Returns
-------
None
"""
rows = int(np.ceil((len(feats))/cols))
if rows==1 and len(feats)<cols:
cols = len(feats)
#print("输入%d个特征,分%d行、%d列绘图" % (len(feats), rows, cols))
if kind == 'hs': #hs:hist and scatter
fig, axes = plt.subplots(nrows=rows*2,ncols=cols,figsize=(cols*5,rows*10))
else:
fig, axes = plt.subplots(nrows=rows,ncols=cols,figsize=(cols*5,rows*5))
if rows==1 and cols==1:
axes = np.array([axes])
axes = axes.reshape(rows,cols) # 当 rows=1 时,axes.shape:(cols,),需要reshape一下
i=0
for f in feats:
#print(int(i/cols),i%cols)
if kind == 'hist':
#frame.hist(f,bins=100,ax=axes[int(i/cols),i%cols])
frame.plot.hist(y=f,bins=100,ax=axes[int(i/cols),i%cols])
elif kind == 'scatter':
frame.plot.scatter(x=f,y='SalePrice',ylim=(0,800000), ax=axes[int(i/cols),i%cols])
elif kind == 'hs':
frame.plot.hist(y=f,bins=100,ax=axes[int(i/cols)*2,i%cols])
frame.plot.scatter(x=f,y='SalePrice',ylim=(0,800000), ax=axes[int(i/cols)*2+1,i%cols])
elif kind == 'box':
frame.plot.box(y=f,ax=axes[int(i/cols),i%cols])
elif kind == 'boxp':
sns.boxplot(x=f,y='SalePrice', data=frame, ax=axes[int(i/cols),i%cols])
i += 1
plt.show()
plotfeats(data_train,feats_continu,kind='scatter',cols=6)
分析上图:
LotFrontage、LotArea、GrLivArea、1stFlrSF、2stFlrSF、GarageArea、BsmtFinSF1、TotalBsmtSF: 这几个面积和距离和售价呈明显正相关趋势
LotFrontage:房子到街道的距离,大多在50-100英尺(15-30米)
LotArea:占地面积(包括房屋、花园、前后院……),均值是10516平方英尺(900+平方米)
GrLivArea:地面以上整体面积
1stFlrSF、2stFlrSF: 第1、 2层建筑面积
GarageArea:车库面积
BsmtFinSF1、BsmtFinSF2、TotalBsmtSF:地下室面积,很多房子还有第2个地下室
YearBuilt、YearRemodAdd、GarageYrBlt:从图中可以看出,建造年限对售价虽正相关,但坡度较小,关联度没有上面几个因素大,早点、晚点售价差不多
plotfeats(data_train,feats_numeric_discrete,kind='scatter',cols=6)
分析上图:
MSSubClass:户型,典型的离散型特征.
OverallQual、OverallCond:房屋材料、新旧度、condition等的整体打分 —— 这个因为是人为打分,有可能存在给售价高的打高分,所以需要持怀疑态度
随机变量的方差的性质:
设X,Y是两个随机变量,E(X),E(Y),D(X),D(Y)分别为各自的期望和方差,则有:
$$D(X+Y)=D(X)+D(Y)+2E(X−E(X))(Y−E(Y)) \tag{$1$} $$特别的,当X,Y相互独立时,有:
$$D(X+Y)=D(X)+D(Y) \tag{$2$} $$对比(1)式和(2)式知,X,Y相互独立时还应该有: $$E(X−E(X))(Y−E(Y))=0 \tag{$3$} $$
这意味着当$E(X−E(X))(Y−E(Y))≠0$时,X与Y不相互独立,而是存在一定关系的。
相关系数和协方差:
我们把$E(X−E(X))(Y−E(Y))$拿出来,单独定义一个概念,即协方差,记为$Cov(X,Y)$,即:
$$Cov(X,Y)=E(X−E(X))(Y−E(Y)) \tag{$4$}$$而
$$ρ_{XY}=\frac{Cov(X,Y)}{\sqrt{DX}\sqrt{DY}} \tag{$5$} $$称为随机变量X,Y的相关系数。
计算每列属性的之间的相关系数,我们关注的是其他属性与SalePrice的相关性,将这些相关系数提取出来之后,按降序排序,可见将相关性最强的几个属性是OverallQual、GrLivArea、GarageCars,恰好我们之前关注过,说明我们平常关于房价的经验在这个数据集上也是符合的。
corrmat = data_train.corr() #计算列与列之间的相关系数,返回相关系数矩阵
f, ax = plt.subplots(figsize=(20, 9))
# 热度图 显示数据值的最大值为0.8,正方形网格
sns.heatmap(corrmat, vmax=0.8, square=True,cmap="YlGnBu")
corr=corrmat['SalePrice'].sort_values(ascending=False)
corr
仔细观察可以发现上面的矩阵中,属性远少于原本的79个,这是因为只有数值型的属性参与的了计算,而离散型的属性并没有,所以需要将这些属性进行处理。
这里仅对CentralAir和Neighborhood进行了转换。
from sklearn import preprocessing
f_names = ['CentralAir', 'Neighborhood']
for x in f_names:
label = preprocessing.LabelEncoder()
data_train[x] = label.fit_transform(data_train[x])
corrmat = data_train.corr()
f, ax = plt.subplots(figsize=(20, 9))
sns.heatmap(corrmat, vmax=0.8, square=True,cmap="YlGnBu")
corr=corrmat['SalePrice'].sort_values(ascending=False)
corr
k = 10 # 关系矩阵中将显示10个特征
cols = corrmat.nlargest(k, 'SalePrice')['SalePrice'].index
cm = np.corrcoef(data_train[cols].values.T)
sns.set(font_scale=1.25)
hm = sns.heatmap(cm, cbar=True, annot=True, \
square=True, fmt='.2f', annot_kws={'size': 10}, yticklabels=cols.values, xticklabels=cols.values,cmap="YlGnBu")
plt.show()
cols
我们在相关性最强的前9个属性中选取特征,这前十个特征分别是:'OverallQual'房屋的整体情况评级、'GrLivArea'地面以上面积,、'GarageCars'车库容量、'GarageArea'车库面积、'TotalBsmtSF'地下室总面积、'1stFlrSF'一层总面积、 'FullBath'地下室全浴室数量,、'TotRmsAbvGrd'地面以上房间数、 'YearBuilt'建造年份。
'GarageCars'车库容量和'GarageArea'车库面积中选取GarageCars作为特征。
'TotalBsmtSF'地下室总面积和'1stFlrSF'一层总面积选取TotalBsmtSF作为特征。
最终考虑的特征如下:
| Variable | Segment | Data Type | Comments |
|---|---|---|---|
| GrLivArea | 1 | 0 | 生活面积 |
| TotRmsAbvGrd | 1 | 0 | 总房间数 |
| FullBath | 1 | 0 | 浴室数量 |
| TotalBsmtSF | 1 | 0 | 地下室总面积 |
| GarageCars | 1 | 0 | 车库 |
| YearBuilt | 0 | 1 | 建造年份 |
| OverallQual | 0 | 1 | 总体评价 |
from sklearn import preprocessing
from sklearn import linear_model, svm, gaussian_process
from sklearn.ensemble import RandomForestRegressor
from sklearn.model_selection import train_test_split
import numpy as np
标准化:原始分数减去平均数然后除以标准差;得到的结果是,对于每个属性/每列来说所有数据都聚集在0附近,方差为1。使得预测结果不会被某些维度过大的特征值而主导
train_data = pd.read_csv('./data/train.csv')
test_data = pd.read_csv('./data/test.csv')
将所有的训练数据和测试数据的79个特征按样本连结,便于处理。
all_features = pd.concat((train_data.iloc[:, 1:-1], test_data.iloc[:, 1:]))
all_features.shape
对连续数值的特征做标准化,然后对这些特征中的空缺值用均值0填充。
numeric_features = all_features.dtypes[all_features.dtypes != 'object'].index
all_features[numeric_features] = all_features[numeric_features].apply( lambda x: (x - x.mean()) / (x.std()))
# 标准化后,每个特征的均值变为0,所以可以直接用0来替换缺失值
all_features[numeric_features] = all_features[numeric_features].fillna(0)
all_features.shape
使用one-hot编码,将离散数值转成指示特征。编码后的属性有331个
all_features = pd.get_dummies(all_features, dummy_na=True)
# dummy_na=True将缺失值也当作合法的特征值并为其创建指示特征
all_features = pd.get_dummies(all_features, dummy_na=True)
all_features.shape
将处理过的数据通过values属性得到NumPy格式的数据,并转成NDArray方便后面的训练。
#划分为训练子集和测试子集
n_train = train_data.shape[0]
train_features = nd.array(all_features[:n_train].values)
test_features = nd.array(all_features[n_train:].values)
train_labels = nd.array(train_data.SalePrice.values).reshape((-1, 1))
X_train,X_test, y_train, y_test = train_test_split(np.array(all_features[:n_train]),np.array(train_data.SalePrice.values).reshape(-1,1), test_size=0.33, random_state=42)
print(X_test.shape)
print(test_features.shape)
print(type(X_test))
print(type(test_features))
test_features = test_features.asnumpy()
print(type(test_features))
由上面结果选择随机森林回归算法
loss = gloss.L2Loss() #L2 loss
def log_rmse( features, labels):
# 将小于1的值设成1,使得取对数时数值更稳定
# clipped_preds = nd.clip(net(features), 1, float('inf'))
# clipped_preds = nd.clip(features, 1, float('inf'))
rmse = nd.sqrt(2 * loss(features.log(), labels.log()).mean()) #❓为什么要根号内要两倍L2loss的平均值
return rmse.asscalar() #大小为1的数组转标量
clfs = {
'svm':svm.SVR(), # 支持向量回归模型
'RandomForestRegressor':RandomForestRegressor(n_estimators=400),#随机森林回归算法
'BayesianRidge':linear_model.BayesianRidge() # 贝叶斯回归模型
}
for clf in clfs:
try:
clfs[clf].fit(X_train, y_train)
y_pred = clfs[clf].predict(X_test)
y_test = nd.array(y_test).reshape((-1, 1))
y_pred = nd.array(y_pred).reshape((-1, 1))
# print(clf + " cost:" + str(np.sum(y_pred-y_test)/len(y_pred)) )
print(clf + " cost:" ,log_rmse(y_pred,y_test))
except Exception as e:
print(clf + " Error:")
print(str(e))
print(X_test.shape)
只选择出相关性大的特征进行随机森林,rmse竟然变大了。
cols = ['OverallQual','GrLivArea', 'GarageCars','TotalBsmtSF', 'FullBath', 'TotRmsAbvGrd', 'YearBuilt']
x = data_train[cols].values
y = data_train['SalePrice'].values
X_train2,X_test2, y_train2, y_test2 = train_test_split(x, y, test_size=0.33, random_state=42)
clf = RandomForestRegressor(n_estimators=400)
clf.fit(X_train2, y_train2)
y_pred = clf.predict(X_test2)
y_test = nd.array(y_test2).reshape((-1, 1))
y_pred = nd.array(y_pred).reshape((-1, 1))
print("RandomForestRegressor" + " cost:" ,log_rmse(y_pred,y_test))
# print(y_pred)
# 保存clf,共下面计算测试集数据使用
# rfr = clf
验证测试集
import pandas as pd
from sklearn.ensemble import RandomForestRegressor
# 之前训练的模型
# rfr=clf
rfr = RandomForestRegressor(n_estimators=400)
rfr.fit(X_train, y_train)
y_te_pred = rfr.predict(test_features)
print(y_te_pred)
print(y_te_pred.shape)
print(x.shape)
prediction = pd.DataFrame(y_te_pred, columns=['SalePrice'])
result = pd.concat([test_data['Id'], prediction], axis=1)
result.to_csv('./data/predictions.csv', index=False)
result.columns
train_data = pd.read_csv('./data/train.csv')
test_data = pd.read_csv('./data/test.csv')
将所有的训练数据和测试数据的79个特征按样本连结,便于处理。
all_features = pd.concat((train_data.iloc[:, 1:-1], test_data.iloc[:, 1:]))
all_features.shape
对连续数值的特征做标准化,然后对这些特征中的空缺值用均值0填充。
numeric_features = all_features.dtypes[all_features.dtypes != 'object'].index
all_features[numeric_features] = all_features[numeric_features].apply( lambda x: (x - x.mean()) / (x.std()))
# 标准化后,每个特征的均值变为0,所以可以直接用0来替换缺失值
all_features[numeric_features] = all_features[numeric_features].fillna(0)
all_features.shape
使用one-hot编码,将离散数值转成指示特征。编码后的属性有331个
# dummy_na=True将缺失值也当作合法的特征值并为其创建指示特征
all_features = pd.get_dummies(all_features, dummy_na=True)
all_features.shape
将处理过的数据通过values属性得到NumPy格式的数据,并转成NDArray方便后面的训练。
n_train = train_data.shape[0]
train_features = nd.array(all_features[:n_train].values)
test_features = nd.array(all_features[n_train:].values)
train_labels = nd.array(train_data.SalePrice.values).reshape((-1, 1))
定义比赛用来评价模型--对数均方根误差Root-Mean-Squared-Error (RMSE):给定预测值 $\hat{y}_1$,…,$\hat{y}_n$ 和对应的真实标签 $y_1$,…,$y_n$ ,它的定义为 $$\sqrt{\frac{1}{n}\sum_{i=1}^n\left(\log(y_i)-\log(\hat y_i)\right)^2}$$
其中使用了L2损失函数,L2-loss的收敛速度要比L1-loss要快得多。
$$L_2=\left|{f(x)-y}\right| ^2$$loss = gloss.L2Loss() #L2 loss
def log_rmse(net, features, labels):
# 将小于1的值设成1,使得取对数时数值更稳定
clipped_preds = nd.clip(net(features), 1, float('inf'))
rmse = nd.sqrt(2 * loss(clipped_preds.log(), labels.log()).mean()) #❓为什么要根号内要两倍L2loss的平均值
return rmse.asscalar() #大小为1的数组转标量
这里选择了基本的线性回归模型。
def get_net():
"""
返回线性回归模型
"""
net = nn.Sequential()
net.add(nn.Dense(1))
net.initialize()
return net
def train(net, train_features, train_labels, test_features, test_labels,
num_epochs, learning_rate, weight_decay, batch_size):
"""
使用线性回归模型去拟合,并使用adam来优化参数。
返回模型在训练集和测试集上的损失值。
net:训练模型
train_features:训练属性
train_labels:训练标签
test_features:测试属性
test_labels:测试标签
num_epochs:优化次数
learning_rate:学习率
weight_decay:权重衰减
batch_size:批量大小
"""
train_ls, test_ls = [], [] #损失值
train_iter = gdata.DataLoader(gdata.ArrayDataset(
train_features, train_labels), batch_size, shuffle=True)
# 这里使用了Adam优化算法
trainer = gluon.Trainer(net.collect_params(), 'adam', {
'learning_rate': learning_rate, 'wd': weight_decay})
for epoch in range(num_epochs):
for X, y in train_iter:
with autograd.record():
l = loss(net(X), y)
l.backward()
trainer.step(batch_size)
train_ls.append(log_rmse(net, train_features, train_labels))
if test_labels is not None:
test_ls.append(log_rmse(net, test_features, test_labels))
return train_ls, test_ls
K 折交叉验证
def get_k_fold_data(k, i, X, y):
"""
返回第i折交叉验证时所需要的训练和验证数据。
"""
assert k > 1 # 如果false,那么raise一个AssertionError
fold_size = X.shape[0] // k # 向下取整
X_train, y_train = None, None
for j in range(k):
idx = slice(j * fold_size, (j + 1) * fold_size) # 切片
X_part, y_part = X[idx, :], y[idx]
if j == i:
X_valid, y_valid = X_part, y_part
elif X_train is None:
X_train, y_train = X_part, y_part
else:
X_train = nd.concat(X_train, X_part, dim=0)
y_train = nd.concat(y_train, y_part, dim=0)
return X_train, y_train, X_valid, y_valid
def k_fold(k, X_train, y_train, num_epochs,
learning_rate, weight_decay, batch_size):
"""
进行K 折交叉验证,用于调节超参数。
k:交叉验证重复K次
X_train:训练集属性
y_train:训练集标签
num_epochs:优化次数
learning_rate:学习率
weight_decay:权重衰减
batch_size:批量大小
"""
train_l_sum, valid_l_sum = 0, 0
for i in range(k):
data = get_k_fold_data(k, i, X_train, y_train)
net = get_net()
train_ls, valid_ls = train(net, *data, num_epochs, learning_rate,
weight_decay, batch_size)
train_l_sum += train_ls[-1]
valid_l_sum += valid_ls[-1]
if i == 0:
semilogy(range(1, num_epochs + 1), train_ls, 'epochs', 'rmse',
range(1, num_epochs + 1), valid_ls,['train', 'valid'])
print('fold %d, train rmse %f, valid rmse %f'% (i, train_ls[-1], valid_ls[-1]))
return train_l_sum / k, valid_l_sum / k
作图相关函数,其中 y 轴使用了对数尺度。
#作图函数semilogy
def semilogy(x_vals, y_vals, x_label, y_label, x2_vals=None, y2_vals=None,
legend=None, figsize=(3.5, 2.5)):
set_figsize(figsize)
# plt.locator_params('y',nbins=10)
plt.xlabel(x_label)
plt.ylabel(y_label)
# plt.ylim(0,1.5)
plt.semilogy(x_vals, y_vals) #对数
if x2_vals and y2_vals:
plt.semilogy(x2_vals, y2_vals, linestyle=':')
plt.legend(legend)
def set_figsize(figsize=(3.5, 2.5)):
"""Set matplotlib figure size."""
use_svg_display()
plt.rcParams['figure.figsize'] = figsize
def use_svg_display():
"""Use svg format to display plot in jupyter"""
display.set_matplotlib_formats('svg')
超参数选择并训练模型
k, num_epochs, lr, weight_decay, batch_size = 5, 100, 5, 0, 64
train_l, valid_l = k_fold(k, train_features, train_labels, num_epochs, lr,
weight_decay, batch_size)
print('%d-fold validation: avg train rmse %f, avg valid rmse %f'% (k, train_l, valid_l))
def train_and_pred(train_features, test_features, train_labels, test_data,
num_epochs, lr, weight_decay, batch_size):
net = get_net()
train_ls, _ = train(net, train_features, train_labels, None, None,
num_epochs, lr, weight_decay, batch_size)
semilogy(range(1, num_epochs + 1), train_ls, 'epochs', 'rmse')
print('train rmse %f' % train_ls[-1])
preds = net(test_features).asnumpy()
test_data['SalePrice'] = pd.Series(preds.reshape(1, -1)[0])
submission = pd.concat([test_data['Id'], test_data['SalePrice']], axis=1)
submission.to_csv('submission.csv', index=False)
train_and_pred(train_features, test_features, train_labels, test_data,
num_epochs, lr, weight_decay, batch_size)